Write a highly optimized custom CUDA kernel for the Softmax activation function, specifically tuned for MXC500 GPU with 1389 GB/s memory bandwidth.

The Softmax function is defined as:

Softmax(x_i) = exp(x_i) / sum(exp(x_j)) for all j

To improve numerical stability, the implementation should use the max-subtraction trick:

Softmax(x_i) = exp(x_i - max(x)) / sum(exp(x_j - max(x)))

Key optimization strategies for MXC500 GPU:
1. **8-element vectorization**: Utilize 8-element vector loads/stores to maximize memory bandwidth utilization
2. **Hierarchical optimization**: Different kernels for different sequence lengths (<=128, 129-2048, >2048)
3. **Warp-level optimization**: For small sequences, use single warp processing
4. **Block-level optimization**: For medium sequences, use block-level reduction
5. **Multi-pass reduction**: For large sequences, use tile-based multi-pass reduction

You should fuse the following steps into a single CUDA kernel:
1. Find the maximum value in each row of the input tensor using vectorized operations
2. Subtract the maximum value from each element in the row
3. Compute the exponential of the result with vectorized operations
4. Sum the exponentials for each row using efficient reduction
5. Divide each exponentiated element by the sum

You are given the following architecture:

import torch
import torch.nn as nn

class Model(nn.Module):
    def __init__(self, dim=-1):
        super(Model, self).__init__()
        self.dim = dim
    
    def forward(self, x):
        return torch.softmax(x, dim=self.dim)

Performance targets for MXC500 GPU:
- Small sequences (<=128): Achieve near-peak warp utilization
- Medium sequences (129-2048): Maximize block-level parallelism
- Large sequences (>2048): Optimize for memory bandwidth utilization
- Overall: Target 2-5x speedup over PyTorch implementation